Strong-field ionization of diatomic molecules and companion atoms: strong-field 
approximation and tunneling theory including nuclear motion 
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We present a detailed comparison of strong-field ionization of diatomic molecules and their com- 
panion atoms with nearly equal ionization potentials. We perform calculations in the length and 
velocity gauge formulations of the molecular strong-field approximation and with the molecular 
tunneling theory, and in both cases we consider effects of nuclear motion. A comparison of our 
results with experimental data shows that the length gauge strong-field approximation gives the 
most reliable predictions. 
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I. INTRODUCTION 

When a molecule is subject to an intense laser field, 
a series of related strong-field processes may occur, in- 
cluding ionization, dissociation and high-harmonic gen- 
eration. As for atoms, it is essential to obtain a detailed 
understanding of the initial, single-ionization process in 
order to describe the subsequent evolution of the system. 
Fully ab initio calculations of the strong-field ionization 
of any molecule more complicated than H 2 are impossi- 
ble in the foreseeable future. Hence, theories of general 
practical use have to rely on theoretical modelling and 
it is the purpose of the present work to investigate the 
accuracy of such models. While the quality of ab ini- 
tio calculations may be checked numerically by studies of 
convergence and by identity of results in different gauges, 
the quality of an approximate model has to be checked 
against a more accurate model or experimental data. 

In the study of strong-field ionization of molecules, 
three models are widely used. These are the velocity 
gauge (VG) molecular strong-field approximation (MO- 
SFA VG) Q, the length gauge (LG) MO-SFA (MO-SFA 
LG) 0, and the molecular tunneling theory The 
latter theory is an extension of the atomic Ammosov- 
Delone-Krainov (ADK) tunneling theory 4] generalized 
to take into account the non-spherical symmetry of the 
molecular system and it is referred to as the MO- ADK 
theory. The MO-SFA VG and MO-SFA LG theories 
are the velocity gauge and length gauge versions of the 
Keldysh-Faisal-Reiss (KFR) atomic theories 0, El ap- 
propriately modified to the molecular case. In short, the 
MO-SFA theories are based on an 5-matrix formulation 
where one considers the transition from a field-free ini- 
tial state to a Volkov final state, i.e., the state of a free 
electron in the laser field. 

To assess the quality of the models, predictions have 
been compared with experimental data (see, e.g., Ref. |3( 
and references therein). In particular, data has been 
studied in detail which compares ionization yields of di- 
atomic molecules with yields of atoms having nearly the 
same electronic binding as the molecules - the so-called 
companion atoms. The ratio of ion signal data of di- 
atomic molecules and their companion atoms 0, 0, 0, 



[Til fl^L Ir3| is illuminating since effects are factored out 
which depend only on the binding energy. Hence a com- 
parison of the ratio of ionization signal of, e.g., N2 and 
its companion Ar atom allows one to study fairly directly 
effects of molecular symmetry and ro- vibrational motion. 
In Ref. Q the status of theory versus experiments was 
discussed: The MO- ADK theory Q gave results in sat- 
isfactory agreement with experiments for molecules with 
suppressed ionization: D 2 :Ar [Tol[T^ |. 2 :Xe HEl and 
without: N 2 :Ar [uj. 

For the F 2 :Ar ratio, the MO- ADK theory predicts sup- 
pression of the molecular signal while the experimental 
data E2 does not. For the F 2 molecule, however, the 
total energy obtained from a Hartree-Fock (HF) calcula- 
tion of the molecular system is higher than the energy of 
the separated atoms [l4|. This means that electron cor- 
relation effects are important for the proper description 
of the binding of F 2 , and, accordingly, the one-electron 
model taken as the starting point for the MO- ADK (and 
the MO-SFA) theory, is not applicable. The models dis- 
cussed above, and in the present paper, can only be 
meaningfully applied on molecular systems where effects 
of electron-electron correlation are small Ej - Methods 
which do take electron-electron correlation into account 
are presently under development [l6l Il7| but these meth- 
ods are not yet applicable to molecules interacting with 
oscillating fields in geometries without cylindrical sym- 
metry. 

The MO-SFA VG was compared with experiments in 
Ref. Q. The theory was shown to predict suppressed 
ionization of 2 :Xe 0, and to predict 'absence' of 
suppressed ionization in N 2 :Ar EH- I n Ref- 0- the 
predictions of the MO-SFA VG were compared with ex- 
perimental above-threshold-ionization spectra of 2 (7T S 
symmetry) and N 2 (o~ g symmetry) and for these systems 
the theory gave the correct qualitative predictions. 

In general, the MO-SFA VG approximation predicts 
suppressed ionization for diatomic molecules with an an- 
tibonding highest occupied molecular orbital (HOMO) 
and it predicts no suppression of the ionization of di- 
atoms with a bonding HOMO. These general predictions 
are based on an analysis of the MO-SFA VG rate in terms 
of linear combinations of atomic orbitals where the an- 
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tibonding [bonding] rate is proportional to sin 2 (q ■ R/2) 
[cos 2 (q ■ -R/2)] and q fl< 1 where q is the momentum 
of the outgoing electron and R the internuclear coordi- 
nate. For example, the MO-SFA VG approximation pre- 
dicts suppressed ionization for O2 with a HOMO of ir g 
symmetry and no suppression of the D2:Ar signal. This 
latter prediction is in contrast with experimental obser- 
vation [13 and MO- ADK theory |3J] , and we shall return 
to this discrepancy below. 

In view of the above discussion at least two questions 
need to be investigated in more detail. Firstly, what is 
the quantitative prediction of the MO-SFA in the D2:Ar 
case? Secondly, is it possible to carry out a calcula- 
tion in a more accurate approximation which sheds some 
light on the discrepancy between the MO- ADK calcula- 
tions of Ref. [3 and the experimental results of Ref. 13] ? 
To address related questions we have recently developed 
the necessary tools for computations in the velocity and 
length gauge formulations of the MO-SFA, and in addi- 
tion we have set up a programme for the evaluation of 
the molecular tunneling theory 0, 0] • Our LG and VG 
versions of the MO-SFA are generalizations of the atomic 
adiabatic theory [ljj. Our formulation of the MO- ADK 
theory follows Ref. 3] , but our method of extracting the 
angular coefficients of relevance for the evaluation of the 
rate is different [2lll5|. In this work, we extend our previ- 
ous analysis, by taking effects of nuclear motion explicitly 
into account and we show that the inclusion of nuclear 
vibrations may lead to a significantly lower rate. For the 
molecules considered in this paper nuclear motion is of 
significant importance for D2 and NO while it leads to 
small effects in the final results for the other molecules 
studied. 

The paper is organized as follows. In Sec. [HI we de- 
scribe the theories and provide a qualitative discussion of 
the effects of nuclear motion. In Sec. IIIII we give some 
calculational details, in Sec. If VI we present the results of 
our calculations, and in Sec. [V] we conclude. 



II. THEORY 



In this section we describe the MO-SFA and MO- ADK 
theories and discuss how to account for effects of nu- 
clear motion. In the MO-SFA one considers a state to 
state transition where the molecular states are generated 
from moving nuclei. Contrary, in the MO- ADK there is 
no detailed specification of the final state in the strong- 
field ionization process, and in order to maintain a qua- 
sistatic tunneling picture one must fix the nuclei at an 
internuclear distance R and let the electron move in the 
independent potential. 



A. Molecular strong-field approximation including 
nuclear motion in the Born-Oppenheimer 
approximation 

In the Coulomb gauge and in the dipole approximation 
the linearly polarized laser field may be described by the 
vector potential A(t) = Aq coswi, where ui is the angular 
frequency. The corresponding electric field is obtained as 
F(t) = —d t A(t), i.e., F(t) = F a smujt. The interaction 
between the field and an TV-electron system is [atomic 
units (h = m e = e = uq = 1) are used throughout] 

y( VG )(t) = f> W . Pj + ^M, (i) 

3=1 

in the velocity gauge and 

N 

vW(t) = '£r j .F(t), (2) 

3=1 

in the length gauge. In either gauge, we express the an- 
gular differential dW/dq and total W rates as sums over 
n-photon absorptions [19j 

= 2tt I A m I 2 q n , (3) 

" n— no 



w = 2?r £ / 1 A i n 1 2 qnd ^ w 

n—no 

where the transition amplitudes corresponding to the in- 
teraction V {c \t), (c = {VG, LG}), 

A^ = ^£{^ f \V^(t)\^)dt, (5) 

involves integration over one period of the field T = — , 
and (||) designates integration over both the electronic 
and nuclear coordinates. Here describes the molecular 
initial state and $>f is the final state describing the state 
of the residual ion and the free electron in the laser field. 
In Eqs. 0-®, n-o is the minimum number of photons 
needed to reach the continuum, and the momentum q n is 
given by energy conservation. In the Born-Oppcnhcimcr 
(BO) approximation q n is determined by Eq. (|10|l below. 

In the SFA we approximate the initial state by a field- 
free molecular state. In the BO approximation this state 
is a product of an electronic state and a ro-vibrational 
state labelled by Ji. The electronic and vibrational 
states are typically the respective ground states. The 
rotational periods of the diatomic molecules are much 
longer than typical experimental pulse durations and 
therefore the rotational degrees of freedom may be ne- 
glected. The total energy of the initial state is 

Ei = Ef(Ro) + E Ui , (6) 
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where E£(Rq) is the electronic eigenenergy at the inter- 
nuclear equilibrium distance Ro and E Vi is the vibra- 
tional eigenenergy of the nuclear Hamiltonian. If we ap- 
proximate the electronic part of the initial wave function 
by the single-determinant HF wave function, the corre- 
sponding initial molecular wave function is 



*i = 



1 



det | il>i(ri)ip 2 (r 2 )...if)N(rN) |n (7) 



where \ Vi (R) is the initial vibrational wave function and 
the ipj's are orthogonal single-electron wave functions. 
The electronic wave function is evaluated at the nuclear 
equilibrium distance Rq since we assume, consistently 
with the BO picture, that it will be a slowly varying 
function of the internuclear distance. We have checked 
that the results are insensitive to this approximation. 

We seek the transition amplitude to a single-electron 
Volkov state and a definite vibrational and electronic 
eigenstate of the molecular ion. The application of a 
Volkov wave in the final state means that the electron- 
ion interaction is neglected. Additionally, we assume that 
the electronic state of the ion is unrelaxed, i.e., only the 



J 



HOMO is affected. The final state is then 
1 



/ 



det | ^iiri)il)2{r2)...^v{rN,t) \r (8) 



x xt{R)e > , 



where ipv is a (27r) -3 / 2 normalized Volkov wave function 
and where the superscripts "+" denote the ionic state. 
The time-averaged energy of the electron in the laser field 
is <? 2 /2 + Up, and the total final-state energy is 



E f =E e / + (R ) + E+ +f + U p , 



(9) 



where U p = Fq/Auj 2 is the ponderomotive potential, and 
where the final state momentum q n is determined by en- 
ergy conservation ntu — Ef — E$, i.e., 

/2{nu - (E e / + (R ) + E+ - Ef(R ) - E Vi ) - U p ). 

(10) 

The transition amplitude of Eq. (JSJ can now be written 
as (c = {VG,LG}) 



4n = \ I (xt f (R)Tpv(r N ,t) | V^it) | 4> N (r N ;R )x„ t (R))^P [% (e £ / + (Ro) + E+ - Ef(R ) - E^ t\dt, (11) 



r 



where the TV-electron matrix element of the one-electron where S Uf:Vi is the Franck-Condon (FC) factor 
operators, Eqs. and was been simplified by the 
Slater-Condon rules [2(j and where Vjif\t) is the tran- 



sition operator for electron V. The integration over nu- 
clear coordinates can be performed immediately to give 



s 1 



ipv( r N, t)V^ c \t)iJj N (r N ;R )dr N 



x e 



X +(R) Xu t (R)dR. (13) 



From Eq. (|12[) . we find the following explicit expressions 
(12) for the amplitudes 



4 L « G) =^4r (- Eb 



[q n + A(t)} \ ^ qn + A ^ expi ( ntJ t + q n . ao sm(ujt) + ^-sm(2ujt) ) I i 1.4) 



A ( q n G) = S Vf>Vi (-E b - M ^Ar(Q„)i J exp i (nut + q n ■ a fan(ut) + ^ sin(2^ t, 



(15) 



In Eqs. (|14fl - i|15|) . a o = Aq/uj is the quiver radius, space wave function of the HOMO, and Eb is the energy 
tPn(q) = (27r)~ 3 / 2 / dre~ iq ' r ipN(r) is the momentum 
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difference between the final and initial state 



Eb — E e f ' + (R ) 



E v - Ef(Ro) — E Vi 



(16) 



When we compare Eqs. (|14|l and 1)15(1 . we see that the 
length gauge formulation accounts for a superimposed 
quiver motion in momentum space of the bound state 
electron via the presence of the A(t) term in the argu- 
ment of the momentum space wave function. Such an 
effect is not present in the velocity gauge amplitude. 

To account for the Coulomb interaction between the 
outgoing electron and the residual molecular ion one has 
typically introduced a Coulomb correction factor. In 
the velocity gauge, this factor is = (k 3 / i^ ) 2Z,on ^ K 

with k = \/2Eb 1] and Z- lon the charge of the residual 

ion, while in the length gauge ^coul = (2k 3 / Fq) 2Z, ° d ^ k 
plf. Both correction factors were derived for the case of 
strong-field ionization of atoms, and hence do not take 
into account the molecular symmetry. In our evaluation 
of rates, we have found that much more precise results 
are obtained in the length gauge with = 1 @> E3 • 

We explain the absense of a Coulomb interaction in the 
length gauge by the fact that the transition to the contin- 
uum occurs at large distances. In this spatial region the 
laser-electron interaction is stronger than the electron-ion 
interaction and the Volkov state is a good approximation 
for the final state. In addition to the Coulomb correction 
factor, we multiply the rates by the number of equiva- 
lent electrons in the HOMO. Finally, to obtain the total 
ionization rate of the molecules, we must sum the con- 
tributions from each vibrational level in the final state. 
In the velocity gauge our result is equivalent to that of 
Ref . • For the noble gas atoms with filled p shells we 
sum the rates from each magnetic sub-state to obtain the 
total rate of ionization. 



1. Qualitative discussion of the effect of nuclear motion 

The transition amplitude of Eq. (|12fl consists of two 
factors, namely the Franck-Condon (FC) factor and an 
electronic matrix element and both factors depend on the 
vibrational levels considered. The rates to each vibra- 
tional level are therefore not just proportional to the FC 
factors. Instead the relative populations in the lower final 
vibrational states are enhanced compared with the dis- 
tribution obtained from the FC factors alone j2^| because 
the electronic matrix element is favoured by the smallest 
energy differences. When including vibrations, the total 
rate summed over all final vibrational states will therefore 
typically be smaller than if the vibrational ground state 
of the ion had been given the weight of unity. This latter 
method is nearly equivalent to fixing the nuclei [compare 
Eqs. (O-O with Eqs. (5)-(6) of Ref. 0]. 

The importance of the inclusion of nuclear vibrations 
will depend on the properties of the neutral molecule 
and the molecular ion. If their potential curves are only 
shifted with respect to each other but otherwise exactly 



identical then the vibrational eigenstates will be identical 
too. The orthogonality of the nuclear wavefunctions then 
assures that only a single FC will be nonzero. We may es- 
timate the importance of nuclear vibrations using molec- 
ular orbital theory and by considering the type of valence 
orbitals. If the valence orbital is nearly non-bonding as, 
e.g., in N2, the bonding properties of the molecular ion 
will be approximately equal to the bonding properties of 
neutral molecule and transitions between the vibrational 
ground states of the molecule and the ion dominate. In 
the case of a bonding HOMO, e.g., as in D 2 , the bond- 
ing of the ion will be weakened and transitions to many 
vibrational states will occur. 



B. Molecular ADK theory 

The tunneling theory Q relies on the assumption 
that at any given instant of time the system will respond 
to the external laser field as if it were a static electric field. 
The rate of ionization in the oscillating field may then be 
determined by time-averaged static rates. Whether this 
approach is reasonable or not depends on the value of 
the Keldysh parameter 7 = y / 2Ebjr |f| with 7 <C 1 in 
the tunneling regime. We will only snow results from the 
MO- ADK theory in the intensity regions corresponding 
to 7 < 1. 

The tunneling rate of diatomic molecules can be de- 
termined once the field-free asymptotic wave function is 
known. In a body-fixed frame, labeled by superscript B 
and a z axis directed along the internuclear axis, this 
function has the asymptotic Coulomb form 

^(r) ~ r z '-/"-V Kr ^C im y im (r), (17) 
1 

where m is the projection of the angular momentum on 
the internuclear axis, and where C\ m are expansion coef- 
ficients. 

From the asymptotic form of Eq. I|17|) , the total ioniza- 
tion rate in a static (DC) field in the positive laboratory- 
fixed Z direction is calculated as in the atomic case 
[H IHil, and the result is @ 



\B(m 



Ml 2 



■f 2\ m '\\m'\\K 2Z ''-/ K - 1 



(18) 



2k 2 



3 \ 2Z ion /K-|m'|-l 



exp 



2 K 
3FoJ ' 



where 



B(m') = 5>l)d m ' 



W)/a (2l + l)(l + \m'\y. 
V 2{l-\m!\)\ 



(19) 



Here d$, (9) is the middle term of Wigner's rotation 
function [25j and 9 is the angle between the field direction 
and the internuclear axis. 
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In a slowly varying field, the ionization rate is found 
by averaging the DC rate over an optical cycle 0] 

w ^ fm w+M) + w- t&i {F Q ) 

V 7TK 3 2 

where W^ at (Fo) are the DC rates for the positive and 
negative field directions with respect to the Z direc- 
tion. When the field points in the negative Z direction, 
Eqs. I |18|l and l|19|l must be modified by the substitution 
Cim -> {-l) l C lm 0. We see from Eqs. (JTSJl-© that 
one only needs to know the C; m coefficients in order to be 
able to evaluate the tunneling formula analytically. See 
also Ref. |l5| for a generalization of the tunneling theory 
to molecules with more than two nuclei. 

The tunneling theory can be extended to include effects 
of vibrations. In the quasistatic picture all potentials seen 
by the active electron should be regarded as being static. 
This means that the nuclei are fixed and the height of the 
tunneling barrier will depend on the internuclear distance 
chosen. We then calculate the rate for each value of R 
and weight these i?-dependent rates by the R probability 
distribution obtained from the nuclear wave function Q . 
As an approximate binding energy we take the energy 
difference between the potential curves of the ion and 
neutral molecule as calculated by HF theory. The vibra- 
tional wave functions are constructed from the harmonic 
approximation to experimental potential curves prj . As 
demonstrated in Sec. IIVI the effect of this i?-dependent 
weighting is quite small. Note in passing that recently 
vibrational distributions were measured and calculated 
with tunneling theory [27j . 



III. CALCULATIONAL DETAILS 



Table I: The molecular and atomic properties necessary for 
the evaluation of the present MO-SFA LG and the MO-ADK 
theory. I p is the experimental adiabatic ionization poten- 
tial and Ro is the equilibrium distance |26j|. Furthermore we 
give the angular coefficients C; m from our Hartree-Fock based 
calculation together with the coefficients from Ref. We 
have chosen the origin at the geometrical midpoint. Franck- 
Condon factors and vibrational energies can be found in the 
references indicated after each molecular species. 
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Ro 


Com 


Clm 


C2m 




C*4m 






(eV) 


(A) 














[32] 


15.47 


0.742 


2.44 
2.51 




0.14 

0.06 




0.00 
0.00 [3] 


N 2 K) 


[33] 


15.58 


1.098 


3.46 
2.02 




1.64 
0.78 




0.12 
0.04 [3] 


o 2 K) 


[34] 


12.03 


1.208 






1.04 
0.62 




0.07 
0.03 [3] 


S 2 {if a ) 


[35] a 


9.36 


1.889 






1.46 
0.81 




0.24 
0.07 [3] 


CO (a) 


[34] 


14.01 


1.128 


-3.93 


2.79 


-1.59 


0.31 


-0.09 










1.43 


0.76 


0.28 


0.02 


[3J 


NO (tt) 


m 


9.26 


1.151 




-0.25 
0.22 


0.82 
0.41 


-0.06 
0.01 


0.04 

m 


SO (tt) 


[36] a 


10.29 


1.481 




1.09 
0.41 


-1.25 
-0.31 


0.34 
0.01 


-0.12 

[3J 


Ar (p) 




15.76 






2.51 
2.44 






m 


Kr (p) 




14.00 






2.59 
2.49 






m 


Xe (p) 




12.13 






2.72 
2.57 









"Based on photoelectron spectrum 



In order to determine the angular coefficients C/ m we 
evaluate the ground state wave function in the HF ap- 
proximation. We seek an accurate description of the 
orbitals at large distances with the correct exponen- 
tial behaviour. For this purpose we find in general the 
usual expansion in an atomic Gaussian basis to be inad- 
equate. Instead we solve the HF equations fully numer- 
ically for the diatomic molecules [28;] and for the atoms 
|29l l30| . After having obtained the ground state orbitals 
we project the highest occupied orbital on the spherical 
harmonics and match the resulting radial functions to 
the form C/ m r Zion / K_1 e~ Kr ' treating the angular coeffi- 
cients as fitting parameters. We give the C; m 's obtained 
by this procedure in Table [I] Since orbitals from the 
HF calculation are optimal within the independent par- 
ticle approximation, the HF C/ m coefficients should be 
more accurate than the multiple scattering coefficients 
reported in Ref. 0. 

The knowledge of the coefficents C; m is sufficient to 
evaluate the MO-SFA LG rate accurately 0]. For the 
evaluation of the MO-SFA VG transition amplitudes, we 
make a numeric Fourier transform of the HOMO. 



If we use the experimental ionization potential in the 
calculation of k, the HOMO is not guaranteed to follow 
the correct asymptotic form of Eq. j!7jl. The asymp- 
totic form will be similar but with k hf = y2jeHOMo] 
substituted for k with Ehomo the eigenvalue of the one- 
electron HF (Roothaan) equation. For the highly corre- 
lated F~ ion, for example, the two different values of k 
differ by 20% and the application of a wrong asymptotic 
wave function will introduce a large error when calculat- 
ing the MO-SFA LG ionization rates j^J. For the sys- 
tems considered here the largest difference between the 
experimental k and k hf is 10% for NO and less than 5% 
for the remaining systems. Thus we conclude that corre- 
lation effects do not affect the outer electron significantly 
for the systems considered here. We believe that the 
smallness of the differences between k and k hf justifies 
the application of the single-active electron models. 

With the ionization rates at hand, we calculate the 
ion signals by integrating over the temporal and spatial 
profile of the pulse [TH l3?| . In all experiments discussed 



6 



o 
55 

75 
c 

D) 
W 



4.0 
3.0 
2.0 
1.0 

4.0 
3.0 
2.0 
1.0 
0.0 



(a) 



/ \ 



MO-SFALG 

MO-SFAVG 

MO-ADK 

Guoetal. (1998) □ 



□ □ ° n ° □ D n ° □ n QHa nOQ ° nn H a "n 5 n □ □ □ g n 



(b) 



DeWittet al. (2001) 



§§5 



10' 



10 



15 



Intensity (Wcm ) 



Figure 1: Intensity dependent ratios between the yields of N2 
and Ar ions. In both panels the laser wavelength is 800 nm. 
In panel (a), the pulse duration is 30 fs and the experimental 
data is from Ref. [Tl|. In panel (b), the pulse duration is 
100 fs and the experimental data is from Ref. |T3| . The thin 
lines are calculations with fixed nuclei and the thick lines are 
the corresponding calculations including vibrations. 



here, the ions were collected through a small pinhole near 
the beam waist and the spatial integration should be re- 
stricted to this region. As the signals saturate, the signal 
ratios will approach unity for such a setup. We calculate 
the rates and signals for each molecular orientation and 
average the signals over orientations in order to simulate 
randomly oriented ensembles. 



IV. RESULTS AND DISCUSSION 



N 2 :Ar 



Several experiments have been performed on N2 and 
its companion atom Ar [tj I12I Hgj . and there are 
significant disagreements between the results. The yield 
ratios were estimated to be N 2 :Ar sa 0.7 11], 0.2 [t|, 1 
[il| and 1.7 0, respectively. Previous theoretical calcu- 
lations with the MO-SFA VG model predicted the N 2 :Ar 
ratio to be above unity [l| while MO-ADK calculations 
predicted the ratio to be 0.4 — 0.6. In Ref. .3] it was 
shown that the experimental differences could not be ex- 
plained by differences in pulse lengths and the possibility 
of dynamical alignment. 

In Fig.^we show our calculations for the N 2 :Ar ratio. 
First, we note that the results from each model are nearly 
independent of pulse duration. In the length gauge we 
predict the ion yields of N 2 and Ar to be nearly identi- 




Intensity (Wcm ) 



Intensity (Wcm" ) 



Figure 2: The ionization rates of N2 and Ar (dotted) at inten- 
sities around the first channel closings obtained from (a) MO- 
SFA LG and (b) MO-SFA VG. For N 2 , we give the rates at 
parallel (solid) and at perpendicular geometry (long-dashed). 
The wavelength is 800 nm and the intensities at which the 
channels close for N2 are indicated by the arrows. Ar has 
a slightly higher ionization potential and the corresponding 
intensities lie at bit lower. 



cal. These findings are in reasonable agreement with the 
experimental data from Refs. [ill IT2] ] . The differences 
between our MO-ADK results and the results of Ref. Q 
lie in the values of the angular C; TO cocfficents. We find 
the ratio to be slightly lower than unity in agreement 
with Ref. [H [Fig. [2(a)]. The velocity gauge MO-SFA 
predicts the ion yield of N 2 to be somewhat larger than 
for Ar in accordance with Ref. pj. We find an oscillat- 
ing behaviour of the signal ratio in the velocity gauge 
which is not supported by the experimental data. The 
origin of this artifact lies in the concept of channel clos- 
ing: At low intensities continuum-continuum transitions 
are very weak. The rate of ionization is then typically 
dominated by the rate which originates from the low- 
est number of photon absorptions, no. As the intensity 
increases the ponderomotive potential rises and eventu- 
ally leads to the closing of the no channel. The effect of 
channel closing on the total ionization rate depends on 
the relative importance of the no process compared with 
the excess-photon processes. As the intensity increases 
and a channel closing approaches, the momentum corre- 
sponding to no photon absorptions will become small. If 
the HOMO does not contain an I = component, then 
i/jn — > as q — > 0, and for such systems the MO-SFA 
VG amplitude, Eq. I|15|) . will be suppressed. Near the 
channel closings the contribution from the no channel is 
thus small and the total rate will not be significantly af- 
fected when this channel closes. If on the other hand the 
HOMO contains an I = component, the low momenta 
will generally be favoured and the no channel gives a large 
contribution of the total rate. In this case, the closing of 
a channel will therefore lead to an abrupt decrease of the 
total rate. These differences were previously mentioned 
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Figure 3: Intensity dependent ratios between the yields of D2 
and Ar ions. In both panels the laser wavelength is 800 nm. 
The thin lines are calculations with fixed nuclei and the thick 
lines are the corresponding calculations including vibrations. 
In panel (a), the pul se duration is 100 fs and the experimental 
data is from Ref. |l3|. In panel (b), the pulse duration is 200 fs 
and the experimental data is from Ref. [To|. 



in Ref. [22j- I n the length gauge, the no channel will be 
important regardless of the type of orbital due to the pre- 
cense of A(t) in the argument of the momentum space 
wave function [see Eq. jH} ] and we always find abrupt de- 
creases of the rate across a channel closing. We show the 
various effects described above in Fig. [21 The HOMO of 
N2 is a <jg (I = component) orbital and correspondingly 
we find a local minimum in the rate around the channel 
closings. The highest occupied atomic orbitals of Ar are 
the degenerate p (I = 1) orbitals, and we find local min- 
ima in the length gauge rate while the velocity gauge rate 
increases smoothly near the channel closings. The differ- 
ences between the parallel and perpendicalar geometries 
of N2 in the two gauges was explained in Ref. ■ 



B. D 2 :Ar 

Two experiments have been reported on D2 and Ar 
[ToL and both showed a suppression of the D2 sig- 
nal compared with the Ar signal. As mentioned in the 
introduction the MO-ADK calculations reproduced this 
result pj while the MO-SFA VG model predicts absence 
of suppression because the HOMO of D2 is a bonding a g 
orbital. 

In Fig. |21 we show our calculations and the experimen- 
tal sets of data. We do indeed find absence of suppressed 
ionization when we use the MO-SFA VG model. On the 
other hand the MO-SFA LG and MO-ADK models both 
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Figure 4: Intensity dependent ratios between the yields of O2 
and Xe ions. In all panels the laser wavelength is 800 nm. 
The thin lines are calculations with fixed nuclei and the thick 
lines are the corresponding calculations including vibrations. 
In panel (a), the pulse duration is 30 fs and the experimental 
data is from Ref. [Tl|. In panel (b), the pulse duration is 
100 fs and the experimental data is from Ref. |T^|. In panel 
(c) the pulse duration is 220 fs and the experimental data is 
from Ref. 0. 



correctly predict suppression. At high intensities the ra- 
tio will approach unity as both signals saturate. The 
MO-SFA LG model has some problems in predicting the 
correct intensity at which this saturation occurs, but we 
note that the experiments do not agree on the saturation 
intensity either. The longer pulse length should be equiv- 
alent to a lower saturation intensity but this is clearly 
not the case when comparing the experimental data in 
Figs. E (a) and (b). 

We see that the inclusion of nuclear vibrations reduce 
the MO-SFA ratios by a significant factor compared with 
the fixed nuclei calculations. The origin of this effect was 
discussed in Sec. Ill A II The point is simply that the 
inclusion of vibrational motion will reduce the rate when 
the molecule and the molecular ion have different bonding 
properties and transitions to many different vibrational 
states occur. 
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Oo:Xe 



Suppressed ionization of C>2:Xe has been observed re- 
peatedly H |H E2. Theoretically, this was explained 
in the MO-SFA VG by the antibonding character of the 
■Kg HOMO of O2 Q. Another explanation within the 
MO-ADK model was given in Ref. 0] where the inter- 
pretation was based on the asymptotic charge density of 
the TT g HOMO. At some molecular orientations the elec- 
tronic density will be preferentially perpendicular to the 
polarization axis and the rate of ionization will then be 
very small. Finally in Ref. [39l | it was proposed that nu- 
clear dynamics could be responsible for the suppression. 
In Fig. 01 we show the experimental data together with 
our calculations. Clearly, all experiments and calcula- 
tions show suppressed ionization of O2. Before satura- 
tion effects become important the theoretical MO-SFA 
VG predicts strongest suppression with the ratio below 
0.01 while the ratio in MO-SFA LG is » 0.03. At low 
intensities, the experimental ratios are scattered between 
0.02 — 0.20 which again makes a quantitative comparison 
with theories difficult. In the tunneling regime, beyond 
the intensity of 10 14 Wcm~ 2 , the MO-ADK model pre- 
dicts that saturation effects are already important and 
thereby the degree of suppression is masked. Experimen- 
tally, saturation sets in at much higher intensities around 
2 x 10 14 WcitT 2 , in good agreement with the MO-SFA 
theories. 
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Figure 5: Intensity dependent ratios between the yields of 
CO and Kr ions. In panel (a) the laser wavelength is 800 nm 
and the pulse duration is 100 fs. In panel (b) the laser wave- 
length is 1365 nm and the pulse duration is 80 fs. In both 
panels the experimental data is from Ref. [TsT] . The thin lines 
are calculations with fixed nuclei and the thick lines are the 
corresponding calculations including vibrations. 



D. CO:Kr 

The ion signal ratio of CO:Kr was measured in Ref. 13] 
and found to be around one half. In Figs. El (a) and 
(b) we show the results at 800 nm and 1365 nm, respec- 
tively. Under the experminental conditions at 800 nm, 
the Keldysh parameter exceeds unity. Thus we cannot 
rely on the MO-ADK theory and in Fig. El (a) we show 
only MO-SFA calculations. Both theories predict the ra- 
tio to be slightly larger than one. In the MO-SFA VG 
model, we find that the ratio depends on the intensity in 
a similar way as the N2:Ar ratio, and since the electronic 
structures of CO and Kr are nearly identical to N2 and 
Ar, respectively, we can explain the oscillating behaviour 
by the same channel-closing argument (see Sec . II V AT) . 

At 1365 nm and in the intensity range of Fig. El (b), 
we would expect to be in the tunneling regime. De- 
spite the fact that the tunneling theory should be ap- 
plicable, the MO-ADK model predicts the ratio to be 
somewhat too large. Finally we note that the result of 
our MO-ADK calculation is approximately an order of 
magnitude larger than the result from Ref. 3] . The only 
differences between these two calculations are the coef- 
ficients Ci lfn —o and the time-averaging of the static field 
rate, our Eq. (|20|l as compared with Eq. (10) of Ref. Q. 



E. Molecules without companion atoms 

In Ref. ^3] the ion signal ratios of the pairs S2:Xe, 
NO:Xe, and SO:Xe were measured. Common to these 
pairs is that the ionization potential of Xe is somewhat 
higher than the ionization potentials of the molecules. 
The ratios measured are thus the results of both the 
structural differences (orientation, electronic wave func- 
tions) and the difference between the binding ener- 
gies. The experiments were performed at the wavelength 
800 nm with 7 > 1, Fig. El and at 1365 nm with 7 < 1, 
Fig.0 

The electronic structure of S2 is similar to O2 and 
we should therefore expect the ionization of S2 to be 
suppressed too. Our MO-SFA calculations, shown in 
Fig. El (a), predict the S2:Xe ratio to be much higher 
than measured. Experimentally, the ratio was deter- 
mined to be around unity. This indicates that ioniza- 
tion of S2 would indeed be suppressed if compared with 
a hypothetical companion ion. In Figs.El(b)-(c) we show 
the calculations for SO:Xe and NO:Xe and find a similar 
disagreement between theory and ex per iment for these 
pairs. The two previous attempts [3J, [ljj to explain the 
ratios in this intensity regime overestimated the ratios by 
3 — 5 orders of magnitude. Both calculations were based 
on tunneling models - in Ref. ^| a purely atomic ADK 
model was used and in Ref. Q the MO-ADK model was 
applied. Our present MO-SFA calculations including nu- 
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Figure 6: Intensity dependent ratios between the ion yields 
of (a) S 2 :Xe, (b) SO:Xe and (c) NO:Xe. In all panels the 
laser wavelength is 800 nm, the pulse duration is 100 fs and 
the experimental data is from Ref. 13]. The thin lines are 
calculations with fixed nuclei and the thick lines are the cor- 
responding calculations including vibrations. 



Figure 7: Intensity dependent ratios between the ion yields 
of (a) S 2 :Xe, (b) SO:Xe and (c) NO:Xe. In all panels the 
laser wavelength is 1365 nm, the pulse duration is 80 fs and 
the experimental data is from Ref. IT?. The thin lines are 
calculations with fixed nuclei and the thick lines are the cor- 
responding calculations including vibrations. 



clear motion are significantly closer to the experimental 
data but the agreement is still poor. 

In Fig.0we show calculations and experimental data at 
the wavelength of 1365 nm, i.e., in the tunneling regime. 
The experiments were performed at such a high intensity 
that saturation effects become important in all our cal- 
culations and correspondingly, the ratios are all around 
unity. These predictions are in agreement with the exper- 
iment for S 2 :Xe, Fig. □(a), and NO:Xe, Fig. 0(c). From 
Fig. |7|(b) we see that the SO:Xe ratio is much lower than 
unity. This is quite remarkable since all three ratios in 
Fig. [5] were nearly identical at 800 nm but very different 
at 1365 nm - such wavelength dependencies are of course 
impossible to predict by quasistatic tunneling theories. 



V. CONCLUSIONS 

We have made a detailed study of strong-field ioniza- 
tion of diatomic molecules and their companion atoms. 
In particular we have investigated whether the approxi- 
mate models: MO-SFA LG, MO-SFA VG and MO-ADK 
are able to correctly predict the presence or absence of 
suppressed ionization in certain molecules. Furthermore 



we have considered how nuclear vibrations can be taken 
into account within the adiabatic theory and the single- 
active electron approximation. 

All models did correctly predict absence of suppressed 
ionization for N 2 :Ar and presence of suppressed ioniza- 
tion for 02:Xe. Quantitative comparisons were made 
difficult due to disagreements between different experi- 
ments, and further experiments would therefore be highly 
desirable. 

We find a rather good overall agreement between 
the MO-SFA LG theory and experiments on diatomic 
molecules and their companion atoms, and we believe 
that the MO-SFA LG model accounts quite well for 
the structural dependency on the strong-field ionization 
rates. In general we found the length gauge MO-SFA to 
be in better agreement with experiments than the veloc- 
ity gauge MO-SFA - in particular in the case of D 2 :Ar 
where the MO-SFA VG predicted no suppression of ion- 
ization in contradication with the experimental findings 
and the other theories. 

The MO-ADK predictions were mostly in accordance 
with previously published results using this model 0. 
For the cases N 2 Ar, D 2 Ar and 2 :Xe these predictions 
are also in qualitative agreement with experimental data. 
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The most significant differences between the MO-ADK 
and the MO-SFA LG calculations are the much lower ra- 
tios predicted by the latter theory in the cases 02:Xe 
and CO:Kr. In both cases the results of the MO-SFA LG 
model agreed better with the experiments. This is inter- 
esting in view of recent experiments |4f| which reported 
suppression of strong-field ionization for transition metal 
atoms relative to expectations, but only compared with 
tunneling theory, and not SFA. For the future work it 
would be interesting to check the SFA LG model on the 
transition metals to investigate if the suppression is due 
to a general breakdown of the single-active electron mod- 
els or if it is due only to a failure of the ADK model. 

All the models considered here fail to predict the cor- 



rect ratio between molecules and atoms with different 
ionization potentials. This indicates that the depen- 
dency on the electronic binding energy is not correctly 
accounted for and points to the need for the development 
of improved models, including electron-electron correla- 
tions and exact final states. 
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